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In electronic structure calculations the optimized effective potential (OEP) is a method that treats 
exchange interactions exactly using a local potential within density-functional theory (DFT). We 
present a method using density functional perturbation theory combined with the Hylleraas varia- 
tional method for finding the OEP by direct minimization which avoids any sum over unoccupied 
states. The method has been implemented within the plane-wave, pseudopotential formalism. Band 
structures for zinc blende semiconductors Si, Ge C, GaAs, CdTe and ZnSe, wurtzite semiconductors 
InN, GaN and ZnO and the rocksalt insulators CaO and NaCl have been calculated using the OEP 
and compared to calculations using the local density approximation (LDA), a selection of gener- 
alized gradient approximations (GGAs) and Hartree-Fock (HF) functionals. The band gaps found 
with the OEP improve on the calculated results for the LDA, GGAs or HF, with calculated values 
of 1.16eV for Si, 3.32eV for GaN and 3.48eV for ZnO. The OEP energies of semi-core d-states are 
also greatly improved compared to LDA. 



I. INTRODUCTION 

Density functional theory (DFT) is one of the great 
successes in tackling the many- electron problerrP^. 
Density-dependent exchange-correlation potentials, such 
as the local density approximation (LDA) or generalized 
gradient approximations (GGAs), are useful and sur- 
prisingly accurate for many properties of materials de- 
spite the known deficiencies of these functionals. Gener- 
ally, the LDA slightly underestimates lattice constants, 
while GGAs slightly overestimate them and both under- 
estimate single-particle band gaps in solidsP^ The ex- 
change energy can be calculated exactly using the non- 
local Hartree-Fock (HF) formalism, rather than being ap- 
proximated as part of the LDA or GGA. However pure 
HF grossly overestimates band gaps in solids and is con- 
siderably more expensive computationally because of the 
fully non-local nature of the exchange operator. 

The central theorem of density functional theory^ 
demonstrates that a local potential V(r) is sufficient to 
completely define a many-electron system in the ground 
state. A local potential that includes the effect of the 
exact exchange interaction was first proposed by Sharp 
and HortorP and solved nearly a quarter of century later 
by Talman and ShadwickP This potential is known as 
the optimized effective potential (OEP), because a local 
potential is created to optimally represent the non-local 
HF exchange potential. The total energy of the OEP 
system is variational with the potential in the same way 
the LDA or GGA total energy is variational in the elec- 
tronic density. This allows the construction of equations 
to describe the potential within a Kohn-Sham DFT for- 
malism. The resultant single particle excitation energies 
agree much better with experiment than those calculated 
using the LDA, GGAs and HFP 

The improved OEP description of excited states is a 
consequence of the freedom from self-interaction error 



of both occupied and unoccupied Kohn-Sham stated^. 
Within the LDA and GGAs the potential experienced by 
all states includes some degree of self-interactiorPl In 
HF the exchange term cancels the self-interaction error 
in the occupied states, but no such cancellation is present 
for the unoccupied states. This results in too small an 
electronic band gap using the LDA and GGAs and too 
large in HfE^HSI. 

There is a distinction between the fundamental band 
gap of a material and the corresponding Kohn-Sham gap. 
The Kohn-Sham gap (E^ s ) is defined as the difference 
between the eigenvalues of the highest occupied orbital 
and the lowest unoccupied orbital 

Ef S = e-N+i - £jVj (1) 

where ti is the eigenvalue of the ith orbital and N is the 
total number of electrons in the system. Whereas the 
fundamental band gap (E g ) is denned as the difference 
in the ionization potential and the electron affinity 

Eg = I — A = E(N + 1) — 2E{N) + E(N — 1), (2) 

where / is the ionization potential, A is the electron affin- 
ity and E(N) is the total energy of a N electron system. 
The two are related by 

E g = Ef s + A xc , (3) 

where A xc is the derivative discontinuity in XC energy 
with respect to particle number j- 6 ! 7 ! 15 ! This constant is 
the energy associated with a change in the exchange- 
correlation potential when a infinitesimal charge is added 
to a N electron system. While a comparison between 
the Kohn-Sham and fundamental band gaps would be 
instructive, the value of the derivative discontinuity is 
unknown for most real mate rials. (See however Godby et 
aP and Gorling et alW^) 
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In finite systems, the exchange potential should decay 
as -1/r in the long-range limit for all states. In the OEP 
this decay occurs for all states irrespective of occupancy. 
In HF the non-local potential for occupied states behaves 
correctly but for the empty states decays exponentially 
and in the LDA/GGA the potential for all states de- 
cays exponentially^. Furthermore, the orbitals in the 
OEP each correctly decay with an individual exponent, 
whereas i n HF a ll occupied orbitals decay with the same 
exponent.™^ 

The OEP is potentially as versatile as LDA/GGA 
and HF methods. It could be used instead of HF as a 
foundation for orbital-dependent potentials (so called hy- 
brid functionals) 3 . Furthermore once the local potential 
has been obtained, the calculation of system properties 
is faster than HF based methods, as the computation- 
ally expensive application of the exchange operator is 
avoided. 

While the OEP corre ctly describes many system prop- 
erties (see Grabo et aP^), it is not as widely used as 
the LDA, GGA and HF-hybrid based methods. This 
can be attributed in part to the lack of comparably ef- 
ficient and robust computational schemes to evaluate it. 
In most formulations of the OEP method a sum over 
all excited states is required and truncation to a finite 
sum yields a slowly-convergent series. Furthermore an 
accurate description of these high energy states requires 
large basis sets at high computational expense if a lo- 
cal basis description is used. Despite these deficien- 
cies, some calculations have been done in solids within 
the linear muffin-tin orbital (LMTO) basis setPH2H, the 
FPLAPW methocP as well as in plane-wave pseudopo- 
tential implementations^ 17 * 18 * 23 * ^. The agreement of the 
calculated results in the previous work above with ex- 
periment can be very good (see Ref. [10] ) . However most 
applications have been to semiconductors. For wide-gap 
insulating systems the performance of the OEP can be 
poor, notably for noble-gas solidsP*l See Ref. [TU] for full 
details and further discussion of the OEP method and its 
results. 

The extremely demanding task of calculating the full 
OEP can be reduced by using the mean-field approxima- 
tion of Krieger, Li and Iafrate^ (well-known as the KLI 
approximation) and even in this less precise form alism 
impressively accurate results have been obtained ^ 25 * 28 *. 

The principal aim of this work is to derive and 
demonstrate a variational method of calculating the full 
exchange-only OEP without the need for a sum over all 
unoccupied states of a system. To this end, our method 
for the calculation of the OEP is formulated using density 
functional perturbation theory (DFPT) and the Hyller- 
aas variational principle. As in the case phonon or elec- 
tric field perturbations, only occupied Kohn-Sham or- 
bitals are explicitly included. This method is then ap- 
plied to a range of semiconductors and insulators. 



II. THEORY 

A. The Optimized Effective Potential 

We first summarize the OEP before deriving a varia- 
tional implementation within the density functional per- 
turbation formalism. For the remainder of this article 
we consider only the exchange-only OEP and will sim- 
ply use the term OEP from here on. A starting point 
for finding the local exchange only OEP is the non-local 
Hartree-Fock equation 



f (r) + f cxt (r) + V n (r) + V£(v) #(r) = ef 0f (r), (4) 



where T(r) is the kinetic energy, l / e xt(r) is the external 
potential, Vn(r) is the Hartree potential, Vx(r) is the HF 
exchange potential, i indexes the electronic states, a is 
the spin index and <j>1{ v ) is th e it h orbital with spin a. 
The Kohn-Sham (KS) equation^ that incorporates the 
OEP is given by 



T(r) + V/ oxt (r) + ^(r) (r) = e?# (r), (5) 



where V a is the effective potential fulfilling the role of 
both the Hartree and exchange potentials. The usual 
derivation of the OEP uses a chain rule expansion of 
the derivative of the exchange-correlation energy with re- 
spect to the density^, however here a treatment inspired 
by perturbation theory is used. 

For a potential which differs from the ground 
state (GS) potential by an amount 6V a (r) such that 
V a (r) = V a {r) + 6V a {r) the GS KS orbitals change by 

C(r) #(r) + J dx*V*(x) Jffi) (6) 



where from first order perturbation theory^ the last term 
of the above equation becomes 
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(7) 



The change in the effective potential also gives a first 
order change in the total energy, thus 



AE[V a ] = J dxSVix 



6E[V a ] 
<5F CT (x)' 



(8) 



The OEP energy 15 [V 7 ] is given by the Hartree-Fock 
functional evaluated using the Kohn-Sham OEP orbitals. 
Using the chain rule and the perturbation of the orbitals, 
the derivative of the energy with respect to the potential 
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is given by 
6E[V 



N" 



6V° (x) J 54>f(r) 8V°(x.) 



(9) 



By applying the Hellmann-Feynman theorem (for brevity 
the explicit r dependence of the potentials has been omit- 
ted) 



5E[V° 



T + V cxt + y H + V£ 



C*(r), 



(10) 



and substituting equations ([7|, ( |To| ) into equation ([9]), we 
obtain 



8E[V] 
8V(x) 



N" 



*E E 



C(r)C*(x)C(x) 



i=l o=JV <r +l 
T + Kxt + Vh + l>x 0r(r)+h.c. (11) 



By substitution of the Kohn-Sham equation (eq. [5]) the 
above can be written as^ 



8E[V a ] 
8V a (x) 



-E E 

i=l a=ZV <r +l 

x C(x)C(x) 



h.c, (12) 



which is known as the OEP equation and was first de- 
rived by Sharp and Hortorffl Here it is important to note 
that the requirement to sum over all unoccupied states 
presents a challenge for convergence of the system prop- 
erties. A very large number of unoccupied orbitals must 
be included to ensure adequate convergence and conse- 
quently the calculation is extremely demanding^ 2 . The 
truncation of the infinite sum has also been the subject of 
concerns about the analytic properties of the solutions^. 



B. Applying The Hylleraas Variational Principle 

We present the following formalism to evaluate the 
OEP which avoids an explicit sum over states by borrow- 
ing ideas from density functional perturbation theory^212ZI 
using the Hylleraas variational method. The effec- 
tive potential is obtained by variational minimization of 
£'[y T (r)]. The minimization direction for the potential 
is defined by the functional derivative of the energy with 
respect to the effective potential 



V^r) -> V a (r) - A 



8E[V a ] 
8V°{y) ' 



(13) 



The functional derivative of the energy can alterna- 



tively be written as 

8E[V a ] 
SV° (x) ~ 

where 

oo 

C(x) = - £ 

a=N" + l 



E 

i=l 



<(x) (C(x))' 



h.c, 



C(x)(C|Vh + v^-x>K) 



(14) 



(15) 



and <^(x) is the first order correction to an unper- 
turbed orbital </>^(x). These first order corrections to 
the orb itals w ere named "orbitals shifts" by Kummel and 
Perdew™ 

An alternative method of calculating </>^(x) without 
explicitly including any unoccupied states is to use the 
the Hylleraas variational principle^. We define a second 
order variational functional 



Gtm = m(T+v ext +v° 



efM) 



(16) 



where is the first order correction to the eigenvalues. 
G1[4>1\ is also variational with respect to (f>f under the 
constraint 



(17) 



and using the Hylleraas variational principle it can be 
shown that the orbital shifts which minimize the second- 



order functional also satisfy equation (15). By substitu- 
tion of the projection operator 



P C = I 



N" 

E 

3=1 



mm 



a=N" + l 



(18) 



equation ( 16 ) may be rewritten 



(19) 



N" 

Ei^)(^i)(^h 



v°M 



+ (CKVk + vi v°){i ^ \<%)(<%\M)- 

The exact first order correction to the orbitals (j>? 

that minimizes Gf[4>°] can also be found using the 
Sternheimer-like equation^ 

(T + V ext + V° -e°M) 

+ ( f ~ E + V x V a M) = 0. (20) 

3=1 
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Equation (20 1 can be solved using iterative methods^. 



Replacing the infinite sum over states in equation ( 15 ) 
with an iterative procedure and considering only the oc- 
cupied Kohn-Sham subspace avoids the severe conver- 
gence difficulties of the former method. Instead of deter- 
mining of how many additional states to include in equa- 
tion (15), at a computational cost increasing with their 



number, the convergence problem is transformed to one 
of determining how many iterative cycles are required. 
Consequently a reliable solution may be achieved at con- 
siderably lower expense than summing equation ( 15 ) di- 
rectly. 



By combining equations ( 14 ) and ( 20 1 , the OEP can be 



found variationally while the first order correction to the 
orbitals is obtained by solving the Sternheimer equation 
directly. 



Equations ( 14 ) and (|20 ) are similar to those given by 
Kiimmel and Pefdew^W Our method differs from theirs 
by exploiting the variational character of £?[V CT ] using 
the explicit gradient of equation ( 14 1 to perform a direct 
minimization. Their method involves a self-consistency 
cycle with an update procedure for vxc using a KLI-likc 
expression and involving a division by the density. Our 
direct minimization is conceptually simpler and, because 
of the variational principle is likely to be more numeri- 
cally robust. 



III. IMPLEMENTATION AND CONVERGENCE 

Our procedure for finding the OEP is to solve a double 
nested loop of minimizations. A flow diagram illustrat- 
ing the full procedure is shown in Fig. 1. The inner of 
the two loops represents the solution of the Sternheimer 



equation (equation ( 20 1) to find the first order correction 
to the orbitals. The method is derived from the Baroni 
Green's function solver technique, employing a conjugate 
gradient minimization scheme to find (f>f (r) in a fixed po- 
tential y CT (r)Pt!! 

Using the first order correction to the orbitals and 
the orbitals themselves, the derivative of the total en- 
ergy with respect to the potential can be found from 
Equation ( 14 ) . This gradient is used in the direction-set 



methods to variationally optimize the potential, and its 
magnitude (residual norm) provides an indication of the 
level of convergence of the effective potential. The varia- 
tion of the potential is accomplished via either one of two 
methods; the first is a conjugate gradient scheme^, the 
second method is a modified steepest descent method, 
known as the Barzilai-Borwein (BB) method, which re- 
quires knowledge of the gradient at the current point and 
previous point only^. 

After each step of varying the effective potential, 
Kohn-Sham orbitals in this new, fixed potential are found 
non self-consistently. When the gradient and difference 
in total OEP Kohn-Sham energy per step is smaller than 
a predetermined threshold, the calculation is considered 
to be converged. 
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Iteratively solve Sternheimer equation 






(eq. 20) 






Calculate conjugate gradient search 






direction (7) for (j>°(r) 












Variationally solve eq. 20 for^^r) 






in direction Y 
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Calculate using eq. 14 
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Vary effective potential: 
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current V a (r) 




Calculate desired material properties 
using V a (r) 



FIG. 1: Flow diagram showing procedure for finding the 
OEP. 



5 



The procedure is initialized by first solving the KS 
equation using a local density-dependent exchange- 
correlation potential, in this case the LDA. From this 
self-consistent calculation the initial orbitals and density 
are obtained. Using this density the effective potential is 
then constructed explicitly using the LDA functional for 
the exchange and correlation. Solving the KS equation to 
find the ground state orbitals is accomplished variation- 
ally using a conjugate gradient method.^ The non-local 
HF exchange energy and the Hartree energy are also ini- 
tially calculated using the LDA orbitals. 

The OEP has been implemented w ithin the pseu- 
dopotential plane-wave code, CASTEP^ 21 ^. Tne or- 
bitals, density and potentials are represented on recti- 
linear grids in the usual manner of a plane-wave DFT 
implementation.^ Kohn-Sham orbitals are described on 
reciprocal space grid points G within a sphere bounded 
by the cut-off wavevector, G max and the density and po- 
tentials are nonzero within a sphere of radius 2G max . 
Therefore the grids in both real and reciprocal space used 
to represent the density and potentials have twice the di- 
mensions of the grid used for the orbitals. The Hyller- 
aas minimization scheme is performed explicitly on these 
real space grids by direct variation, so that the effective 
basis used to represent V° (r) is the set of grid points 
{G} : |G| < 2G max . 

Optimized norm-conserving pseudopotentials used in 
this work were generated using the Opium code™ devel- 
oped by Rappe et al. The Hartree- Fock approximation 
was used and the non-analytic behaviour of HF pseu- 
dopotentials was treated using the localization and op- 
timization scheme of Al-Saidi, Walter and Rapped In- 
cluding exchange exactly in the pseudopotentials is re- 
quired to ensure accurate core- valence interaction within 
the current formalism. Incorrect core-valence interaction 
has been shown to have n oticeable effects on calculated 
electronic structures We found that calculations us- 
ing pseudopotentials which include d-states in the core 
typically predict band gaps 20-70% larger than if the d- 
states are treated as valence. Therefore semi core d-states 
were treated as valence for As, Cd, Ga, Ge, In, Se, Te 
and Zn. 

The convergence characteristics of the conjugate- 
gradient and Barzilai-Borwein variational minimizers, 
are compared in Figure 2 which plots the total OEP en- 
ergy for diamond as a function of iteration number. An 
initial rapid decrease in the total energy is followed by a 
long tail of decreasing change in the total energy in both 
cases. The conjugate gradient energy decreases smoothly 
and monotonically towards the converged result, but the 
Barzilai-Borwein energy does not, exhibiting sharp drops 
and an occasional spiked increase (behaviour which has 
been noted previously^ 21 ). After 100 or so steps the BB 
method still converges rapidly while the conjugate gra- 
dient method exhibits a very slow energy decrease. The 
difference in total energy per step is below 2.5/ieV per 
atom after 250 iterations of the BB minimizer and the 
solution is stationary after 300 iterations. For the conju- 
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FIG. 2: Convergence of total Kohn-Sham energy for the 
OEP against number of steps. The solid line is the 
Barzilai-Borwein minimizer and dashed line for the 
conjugate gradient minimizer. Inset: Energy residual on a 
log-linear scale. 



gate gradients method however the energy difference per 
step does not drop below 5/ieV per atom even after 500 
iterations. After 500 iterations the energy difference be- 
tween the two methods is within 8/zeV per atom. In both 
cases we find that the calculated gradient tends to zero as 
the total energy converges, as expected in a variational 
method. The convergence rates of the T point band gap 
of diamond behave in a similar fashion as shown in Fig- 
ure 3. After 250 iterations the gap is converged to within 
2.5/ieV per atom for the BB minimizer. By comparison, 
with the conjugate gradient minimizer the band gap is 
converged to 50/ieV per atom in the same number of it- 
erations. All calculations presented below using the OEP 
method were run for at least 250 iterations to ensure suf- 
ficient accuracy. 



In the calculations that follow, the basis set size (plane- 
wave cut-off energy) and Brillouin-Zone sampling were 
chosen so that total energy differences, evaluated using 
the LDA, were less than 1 meV/atom. The same settings 
were used for the OEP calculations, which resulted in a 
very similar convergence error with cutoff as the LDA 
case. However the OEP calculations were more sensitive 
to Brillouin-Zone sampling error; for example, in the case 
of diamond, 8x8x8 resulted in a 1 meV/atom error using 
the LDA but 3 meV for OEP rising to 9 meV for ZnO. 
Convergence testing was also performed on the FFT grid 
used to represent the potential and the energy error was 
determined to be ss 10/xeV. 
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FIG. 3: Convergence of the F point band gap of diamond 
using the OEP against number of steps. Solid dots are for 
the Barzilai-Borwein minimizer and open dots are for the 
conjugate gradient minimizer. 



IV. CALCULATED ELECTRONIC PROPERTIES 

An example OEP is displayed in figure |4j in this case a 
slice through the primitive unit cell of Si in the diamond 
structure. The potential is smoothly varying outside the 
pseudopotential regions, which are the only regions where 
the potential is positive. The potential is at its deepest 
in the bonding region between the two ions and (outside 
the pseudopotential core region) is always negative but 
tending to zero in the low density regions away from the 
ions. 

The electronic structure of a selection of insulating and 
semiconducting materials were calculated using the LDA, 
some GGAs (PBE 51 , PBESOlP, PW9P 1 and WCP), 
HF and the OEP. The resultant band gaps are shown in 
table I and for the LDA, OEP and Hartree-Fock calcula- 
tions are also plotted in figure [5| The mean absolute error 
is 0.55eV for the OEP, 1.39eV for the LDA and 1.19eV for 
PBE, although these values are skewed by the underesti- 
mated values for the band gap found for diamond, NaCl 
and CaO materials. The corresponding band structures 
for Ge, CdTe, InN, GaN and ZnO are plotted in figures 

mm 

Figure [6] shows the band structure of Ge in the dia- 
mond structure calculated using the LDA and the OEP. 
For Ge the OEP improves the gap from the LDA value of 
0.09eV to 0.86eV which is very close to the experimental 
gap. Similarly for Si the gap is improved from 0.43eV 
using the LDA to 1.16eV using the OEP. For GaAs the 
OEP opens the gap from the LDA value of 0.99eV to 
1.86eV. The OEP gap is greater than the experimental 
gap by 20%; a similar overestimation of the gap by the 
OEP in GaAs was also observed by Gorling et aZ™ Fig- 



ure[7]shows the band structures for CdTe calculated using 
the LDA and the OEP. The OEP gap is 35% larger than 
the experimental gap, opening from 1.46eV for the LDA 
to 2.20eV for the OEP. 

Figure [8] shows the band structures for the wurtzite 
semiconductor InN calculated using the LDA and the 
OEP. The calculated band gap for InN has a small direct 
gap of 0.21eV for the LDA which is opened to a gap of 
1.39eV for the OEP. Similarly the gap in GaAs and CdTe 
is overestimated when compared to the experimental gap. 
The band structures calculated for GaN in the wurtzite 
structure are shown in figure [9j The OEP band gap of 
3.32eV is very close to the experimental gap of 3.39eV 
and greatly improved on the LDA band gap 2.13eV. The 
wurtzite structure ZnO band structures calculated using 
the LDA and the OEP are shown in figure [To} The cal- 
culated band gap is greatly improved from 1.93eV when 
using the LDA to 3.48eV for the OEP which is very close 
to the experimental gap of 3.43eV. 

For ZnSe in the zinc blende structure the OEP cal- 
culated gap of 2.83eV is remarkably close to the exper- 
imental gap and greatly improved over the LDA value 
of 1.88eV. However for diamond the OEP gap of 4.87eV 
underestimates the experimental value by 12% but still 
improves on the LDA value of 3.98eV, the underestima- 
tion by the OEP was again noted by Gorling et alW^ 

For the wide gap rocksalt structure insulators CaO and 
NaCl the OEP band gaps are opened up compared to 
the LDA band gaps, with the LDA calculated band gap 
being 3.93eV for CaO and 4.84eV for NaCl and the OEP 
calculated band gap being 6.09eV for CaO and 6.32eV for 
NaCl. For CaO the gap is underestimated by 14% for the 
OEP and for NaCl the predicted gap is 30% smaller for 
the OEP than the experimental value. The OEP value for 
CaO is very similar to that found by Kotani and AkaP 
and for NaCl is in line with that found by Li et aiP^ 

Table II compares the OEP band gaps obtained here to 
those obtained by others using a variety of basis sets. The 
band gaps calculated here are very similar to previous 
values and confirm the progressive underestimation of the 
band gap for wide gap materials. 

The OEP also improves values for semi-core energy lev- 
els relative to the LDA, as shown in table III. This is to 
be expected given that the self-interaction error present 
in the LDA has the largest effect on the tightly-bound d 
states, and that there is no self-interaction error within 
the OEP. In Ge the 3d electrons occupy states in the 
range -31.9eV to -32.0eV below the valence band max- 
imum with the experimental binding energy being de- 
termined as 29.1eV to 29.6eV. The LDA puts them in 
the range -35.0eV to -35.2eV. The experimental In 4d 
electron binding energy in InN is 17.4eV; the LDA pre- 
dicts between -18.4eV and -18.9eV and the OEP between 
-16.3eV and -17.2eV. 

For GaAs the experimental Ga 3d states lie between 
18.6eV and 19.0eV below the valence band maximum; 
the OEP predicts a range of between -20.0eV and -20.2eV 
and the LDA between -22.9eV and -23.0eV. Experiment 



7 




FIG. 4: The exchange potential (Voep — Vk) as a slice perpendicular to the [111] direction (bond axis) in diamond structured 
Si for the primitive unit cell. Atomic positions are at the origin and quarter of the long diagonal. The colour key is in eV. 
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2.35 


2.23 


2.36 


2.23 


10.38 


3.32 


3.39^ 


ZnO 


1.93 


2.28 


2.09 


2.29 


2.10 


11.51 


3.48 


3.43 s1 


C 


3.98 


4.21 


4.03 


4.24 


4.08 


12.76 


4.87 


5.47 s1 


CaO 


3.93 


4.11 


4.01 


4.12 


4.02 


14.63 


6.09 


7.09 s1 


NaCl 


4.84 


5.34 


5.13 


5.35 


5.19 


13.74 


6.27 


8.97 s1 



TABLE I: Energy gaps (in eV) for some semiconductors and insulators. Calculated values shown for LDA, PBE, PBESOL, 

PW91, WC, Hartree-Fock and the OEP method. 





pp p^H7ll8l62l 


FLAPW™ 


LMTO-ASA 12 ^^ 


KKR-ASA 12 ^^ 


This Work 


Experimental 


Ge 


0.94 


0.89 


1.57, 1.12 


1.03 


0.86 


0.79 


InN 


1.40* 








1.39 


0.93 


Si 


1.14, 1.23 


1.30 


1.93, 1.25 


1.12 


1.16 


1.16 


GaAs 


1.78 


1.74 






1.86 


1.52 


GaN 


2.76*,3.46*,3.49* 








3.32 


3.39 


ZnO 


2.34* 








3.48 


3.43 


C 


5.06, 4.90 


5.20 


5.12, 4.65 


4.58 


4.87 


5.47 


CaO 






6.15 


6.29 


6.09 


7.09 


NaCl 




6.3f 






6.27 


8.97 



TABLE II: Comparison of calculated OEP band gaps with previous work. The basis sets are plane wave pseudopotential 
(PW PP), full potential linearized augmented plane wave (FLAPW), linear muffin tin orbitals with the atomic sphere 
approximation (LMTO-ASA) and Korringa-Kohn-Rostoker method in the atomic-sphere approximation (KKR-ASA). *These 
results are calculated for the zinc-blende structure. fThis result is obtained with the self interaction corrected functional of 
Perdew and Zunger^ 1 rather than the Hartree-Fock exchange functional. 



Material 


State 


LDA 


OEP 


Expt. 


Ge 


Ge 3d 


35.01-35.16 


31.89-32.01 


29.1-29.6^ 


InN 


In 4d 


16.34-17.22 


18.38-18.94 


17 M 


GaAs 


Ga 3d 


22.95-23.01 


20.06-20.21 


18.60-19.04^ 


GaAs 


As 3d 


47.45-47.46 


44.45-44.49 


40. 37-41.0^ 


CdTe 


Cd 4d 


11.89-12.30 


10.16-11.03 


10.1(P 


ZnSe 


Zn 3d 


11.65-11.81 


9.45-9.74 


9.2^ 


ZnSe 


Se 3d 


61.02-61.03 


58.50-58.52 


55.5^ 


GaN 


Ga 3d 


21.23-21.55 


18.30-18.85 


17.94^ 


ZnO 


Zn 3d 


8.88-9.69 


6.65-7.61 





TABLE III: Binding energies (eV) relative to the valence band maxima for semi-core electrons calculated using the LDA and 

OEP. 



Experimental Bandgap (eV) 




FIG. 5: Comparison of calculated and experimental band 
gaps, red squares are HF, blue diamonds are OEP and green 
circles are LDA. (Colour online) 




l r x wk _1 ^v l r x wk 



FIG. 7: Band structure of CdTe (zinc blende structure), 
calculated using the LDA (left) and the OEP (right). 





"W l r x wk ~ l5 w l r x wk 

FIG. 6: Band structure of Ge (diamond structure), 
calculated using the LDA (left) and the OEP (right). 



" 20 r a h k r m l h " 20 r a h k r m l h 

FIG. 8: Band structure of InN (wurtzite structure) using 
the LDA (left) and the OEP (right). 



gives As 3d states at 40.4eV to 4T.0eV below the valence 
band maximum; the OEP puts them at -44.5eV and the 
LDA predicts -47.5eV. In CdTe the Cd 4d electrons lie 
between, -10.2eV and -ll.OeV for the OEP, -11.9eV and 
-12.3eV for the LDA, compared with an experimental 
value of lO.lOeV below the valence band maximum. The 
experimental results also place the Te 5s electrons very 
close to the Cd 4d states with a binding energy of 9.2eV 
below the valence band maximum which the OEP also 
appears to predict. For ZnSe the Zn 3d electrons lie in the 
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FIG. 9: Band structure of GaN (wurtzite structure) using 
the LDA (left) and the OEP (right). 




FIG. 10: Band structure of ZnO (wurtzite structure) using 
the LDA (left) and the OEP (right). 



range -9.4eV to -9.7eV below the valence band maximum 
for the OEP, -11.6eV to -11.8eV for the LDA and 9.2eV 
from experiment. For the Se 3d electrons the LDA gives 
energies of -61.02 compared with -58.5 eV with the OEP 
and 55.5eV from experiment. 

For GaN the LDA gives the Ga d-electrons lying be- 
tween -21.2eV and -21.6eV, the OEP between -18.3eV 
and -18.9eV with an experimental value of 17.9eV below 
the valence band maximum. Using the OEP for ZnO 
gives Zn 3d states in the range -6.6eV to -7.6eV relative 
to the valence band maximum compared with an exper- 



imental value of 7.4eV and a range of -8.9eV to-9.7eV 
with the LDA. As can be seen from the band structures 
for the above materials in the valence, the s and p states 
are almost identical when using the LDA and OEP. It is 
the d states which display the greatest difference. 

The OEP also improves upon the predicted electronic 
structure given by the GGA methods of PBE, PBESOL, 
PW91 and WC which underestimate the gaps of the ma- 
terials investigated here and Hartree-Fock which greatly 
overestimates the gaps. Further work on magnetic metal- 
oxides will be reported in a future publication. 

V. CONCLUSIONS 

A method of obtaining the OEP which treats the lo- 
cal exchange potential exactly without using a sum over 
all unoccupied states has been derived using the Hyller- 
aas variational method and ideas borrowed from density 
functional perturbation theory. This allows for the cal- 
culation of the OEP using a variational minimization 
scheme in real space. The electronic structure of well 
known materials with a wide selection of band gaps have 
been calculated and the band gaps for semiconductors 
are found to be in good agreement with experimental 
values, although for the larger band gap materials, di- 
amond, CaO and NaCl, the calculated band gaps are 
still underestimated by 10-30%. Hartree-Fock pseudopo- 
tentials were found to give more accurate results than 
LDA pseudopotentials. The absence of self-interaction 
error within the OEP is manifest in a better description 
of semi-core d-states compared to the LDA. Their ener- 
gies with respect to the top valence band are much closer 
to experimental spectroscopic measurements than within 
the LDA. 
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